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Abstract 

A  spline-based  method  for  approximating  thin  shell  dynamics  is  presented  here.  While  the 
method  is  developed  in  the  context  of  the  Donnell-Mushtari  thin  shell  equations,  it  can  be 
easily  extended  to  the  Byrne-Fliigge-Lur’ye  equations  or  other  models  for  shells  of  revolution 
as  warranted  by  applications.  The  primary  requirements  for  the  method  include  accuracy, 
flexibility  and  efficiency  in  smart  material  applications.  To  accomplish  this,  the  method  was 
designed  to  be  flexible  with  regard  to  boundary  conditions,  material  nonhomogeneities  due 
to  sensors  and  actuators,  and  inputs  from  smart  material  actuators  such  as  piezoceramic 
patches.  The  accuracy  of  the  method  was  also  of  primary  concern,  both  to  guarantee  full 
resolution  of  structural  dynamics  and  to  facilitate  the  development  of  PDE-based  controllers 
which  ultimately  require  real-time  implementation.  Several  numerical  examples  provide  initial 
evidence  demonstrating  the  efficacy  of  the  method. 
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1  Introduction 


A  spline-based  method  for  numerically  approximating  thin  shell  dynamics  is  presented  in 
this  paper.  This  research  is  motivated  by  the  growing  number  of  applications  involving  shell 
models  for  which  explicit  solutions  are  typically  unavailable.  A  specific  example  involves  the 
development  of  control  strategies  for  reducing  vibrations,  fatigue  and  structure-borne  noise  in 
aircraft  fuselages.  Typically,  some  form  of  cylindrical  shell  model  is  employed  with  the  specific 
equations  and  boundary  conditions  dictated  by  the  application  or  experimental  setup.  When 
numerically  discretizing  the  models  in  such  applications,  the  approximation  method  must  be 
flexible  with  respect  to  boundary  conditions  and  adaptable  with  regard  to  nonhomogeneities 
in  materials  and  geometries.  For  example,  if  piezoceramic  patches  bonded  to  or  embedded  in 
the  shell  are  used  as  actuators  or  sensors,  the  numerical  method  must  be  sufficiently  flexible  so 
as  to  permit  extension  to  models  incorporating  these  components.  Also,  the  primary  motion 
in  many  vibration  and  noise  control  applications  is  bending-dominated  and  the  numerical 
method  must  accurately  approximate  such  dynamics.  Finally,  from  a  control  perspective,  the 
approximation  method  must  adequately  preserve  stability  properties  of  the  physical  system 
as  well  as  be  sufficiently  efficient  so  as  to  permit  real-time  implementation. 

Current  techniques  for  approximating  shell  dynamics  include  modal  expansions  [17,  18, 
26],  Unite  element  discretizations  [2,  3,  4,  11,  13,  14,  15,  20,  21,  22,  25,  28,  30]  and  finite 
difference  approximations  [31,  32].  From  a  theoretical  perspective,  modal  expansions  arise 
naturally  when  separating  variables  in  linear  shell  models.  In  models  for  which  analytic 
expressions  for  the  eigenfunctions  or  modes  can  be  obtained,  this  provides  a  straightforward 
and  often  quite  efficient  method  of  approximation.  The  difficulty,  however,  lies  in  the  fact 
that  analytic  expressions  describing  mode  shapes  are  known  only  for  a  restricted  class  of 
boundary  conditions  and  for  a  limited  set  of  models.  For  general  boundary  conditions  and 
models,  or  systems  involving  coupling  between  a  shell  and  an  adjacent  component  (e.g.,  an 
acoustic  field  or  a  piezoceramic  actuator),  the  modes  must  first  be  numerically  approximated 
or  experimentally  determined  before  expansions  can  be  formed. 

Experimental  determination  of  natural  frequencies  and  modes  is  typically  accomplished 
by  exciting  the  structure  using  an  acoustic  source,  shaker,  magnetic  driver,  et  cetera.  At  each 
natural  frequency,  the  shape  of  the  corresponding  mode  is  determined  and  characterized  for 
use  in  modal  expansions.  As  detailed  in  [32,  33],  however,  difficulties  are  en counted  in  this 
procedure  in  plate  and  shell-like  structures  due  to  the  presence  of  multiple  independent  mode 
shapes  which  can  occur  at  single  experimental  frequencies.  In  such  cases,  the  experimentalist 
must  excite  the  structure  at  various  locations  and  complete  an  orthogonalization  process  to 
obtain  a  complete  modal  basis.  The  determination  of  which  frequencies/modes  to  test  for 
this  behavior  requires  numerical  analysis  or  extreme  care  when  performing  the  experiments. 
Furthermore,  as  discussed  in  [32],  internal  or  air  damping  will  cause  distortion  of  modes  (e.g., 
modal  lines  that  do  not  cross)  which  if  not  accounted  for,  will  degrade  modal  expansions  for 
approximating  the  structural  dynamics. 

In  applications,  modes  determined  for  boundary  conditions  yielding  analytic  expressions 
(e.g.,  simply-supported  end  conditions)  are  occasionally  used  to  approximate  solutions  to  mod¬ 
els  with  more  general  boundary  conditions  or  models  incorporating  additional  components. 
In  some  cases,  the  influence  of  various  boundary  conditions  on  modal  characteristics  has  been 
investigated  [17];  however,  without  any  convergence  analysis  for  such  techniques,  convergence 
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of  the  numerical  approximates  to  the  true  model  solution  cannot  be  guaranteed.  Finally,  it 
is  difficult  to  apply  modal  techniques  in  applications  modeled  by  nonlinear  shell  equations 
or  nonlinear  actuator  contributions  (e.g.,  electrostrictive  actuators  or  piezoceramic  actuators 
embedded  in  constrained  damping  layers). 

Finite  element  approaches,  on  the  other  hand,  are  directly  applicable  for  a  variety  of 
boundary  and  coupling  conditions  and  in  models  including  actuator  contributions  or  nonlinear 
dynamics.  Moreover,  a  large  body  of  software  exists  for  obtaining  finite  element  solutions  for 
various  shell  models  (see  [13,  20,  21,  22]  for  discussions  of  finite  element  methods  and  software 
for  shells  of  revolution). 

A  difficulty  when  developing  and  applying  finite  element  methods  for  shell  applications, 
however,  concerns  the  manifestation  of  various  locking  phenomena.  Shear  locking,  which 
has  also  been  studied  extensively  in  the  context  of  Reissner-Mindlin  plate  models,  is  due  to 
element  incompatibility  when  enforcing  the  Kirchhoff-Love  constraint  of  vanishing  transverse 
shear  strains  as  the  shell  thickness  h  tends  to  zero  [1,  10].  In  shell  applications,  an  even  more 
serious  problem  leading  to  the  failure  of  various  finite  element  methods  is  the  phenomenon  of 
membrane  or  shear-membrane  locking.  This  phenomenon  occurs  when  the  total  deformation 
energy  is  bending  dominated,  and  is  due  to  smoothness  and  asymptotic  constraints  in  the  shell 
model  which  are  not  appropriately  represented  by  the  approximation  method  (see  [3,  4,  15,  25, 
28]).  If  these  constraints  are  not  satisfied  by  approximating  elements,  the  numerical  solution  is 
often  overly  stiff  in  the  sense  that  the  model  exhibits  bending  dynamics  which  the  approximate 
solution  cannot  match.  As  detailed  in  [25],  mesh  sizes  must  be  chosen  significantly  smaller 
than  the  shell  thickness  to  ensure  accurate  approximations  with  high-order  finite  elements  in 
such  bending  dominated  regimes.  These  examples  also  illustrate  that  even  with  such  mesh 
size  restrictions,  low-order  finite  element  methods  often  fail  in  such  regimes.  The  use  of  finite 
elements  which  exhibit  locking  is  detrimental  in  simulations  and  will  seriously  degrade  the 
performance  of  any  controllers  designed  around  such  methods. 

A  second  issue  which  must  be  considered  when  developing  finite  element  methods  for  shells 
concerns  the  piecewise  constant  material  parameters  and  inexact  boundary  conditions  which 
often  arise  in  smart  material  applications.  For  example,  the  use  of  piezoceramic  actuators  leads 
to  piecewise  discontinuities  in  the  density,  stiffness,  Poisson  ratio  and  damping  parameters. 
Consequently,  finite  element  meshes  must  be  aligned  with  the  actuator  boundaries  to  maintain 
accuracy  (the  reader  is  referred  to  [34]  for  a  discussion  of  finite  elements  designed  specifically 
for  piezoelectric  applications).  In  experiments  involving  shells,  slight  energy  loss  must  often  be 
modeled  into  the  boundary  conditions  to  account  for  inexact  boundary  clamps.  Consequently, 
elements  must  be  provided  with  extra  degrees  of  freedom  to  accommodate  these  boundary 
conditions.  Standard  commercial  codes  not  providing  these  capabilities  will  lead  to  potentially 
inaccurate  results  when  applied  to  shells  with  smart  material  actuators  and/or  sensors. 

The  use  of  finite  difference  discretizations  for  approximating  thin  shell  dynamics  is  less 
common  due  to  inherent  difficulties  with  the  high-order  equations  and  boundary  conditions 
arising  in  the  models.  The  reader  is  referred  to  [31]  for  further  discussion  of  finite  difference 
methods  for  shell  applications. 

In  this  paper,  we  consider  a  Galerkin  method  for  discretizing  thin  shell  models  with  linear 
or  cubic  splines  chosen  as  basis  elements  for  approximating  the  longitudinal  and  circumferen¬ 
tial  motion  and  cubic  splines  used  to  approximate  the  transverse  component.  These  choices 
are  motivated  by  smoothness  and  accuracy  criteria,  adaptability  for  a  variety  of  boundary 
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conditions,  and  suitability  when  including  the  contributions  of  actuators  such  as  piezoceramic 
patches  bonded  to  the  shell  [6].  For  the  discussion  here,  the  shells  are  assumed  to  have  lengths 
that  are  relatively  short  in  relation  to  the  radii,  and  the  Donnell- Mushtari  shell  equations  are 
employed.  This  facilitates  the  discussion  while  providing  a  framework  which  is  easily  extended 
to  more  accurate  long  shell  models  through  the  inclusion  of  Byrne-Fliigge-Lur’ye  components. 
In  this  study,  the  formulation  of  the  method  and  numerical  examples  demonstrating  the  ac¬ 
curacy,  efficiency  and  flexibility  of  the  method  are  presented.  Emphasis  in  these  examples 
is  placed  on  demonstrating  that  when  hx  denotes  the  axial  mesh  size,  the  expected  0{h^) 
and  0{h‘^)  accuracy  of  the  method  is  maintained  when  approximating  with  linear  and  cubic 
splines,  respectively.  Convergence  analysis  and  further  analysis  of  the  method  with  regards 
to  membrane  and  shear  locking  will  appear  in  a  future  companion  paper. 

In  Section  2,  the  strong  and  weak  forms  of  the  modeling  equations  for  a  thin  shell  with 
surface-mounted  piezoceramic  patches  are  summarized.  Modal  analysis  for  the  special  case 
of  an  undamped  shell  with  constant  coefficients  and  simply-supported  boundary  conditions  is 
discussed  in  Section  3.  This  provides  a  framework  for  testing  the  convergence  of  the  approxi¬ 
mate  mass  and  stiffness  components  in  the  system.  The  approximation  method  and  resulting 
finite  dimensional  matrix  system  are  detailed  in  Section  4.  Finally,  examples  illustrating  the 
method  are  presented  in  Section  5.  These  include  modal  approximations  for  various  bound¬ 
ary  conditions  as  well  as  results  demonstrating  convergence  rates  when  approximating  shell 
dynamics  generated  by  various  external  inputs.  To  illustrate  the  shell  dynamics  generated  by 
the  piezoceramic  patches,  a  modeled  voltage  spike  to  the  patches  is  used  as  input  in  the  final 
example.  The  resulting  frequencies  are  then  compared  with  those  obtained  in  the  eigenvalue 
(modal)  analysis  to  demonstrate  the  consistency  of  the  method. 


2  Donnell-Mushtari  Equations 


X 


Figure  1.  Thin  cylindrical  shell  with  surface  mounted  piezoceramic  patches. 
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We  consider  here  a  shell  of  length  thickness  h  and  radius  R  having  mass  density  p, 
Young’s  modulus  E,  Poisson  ratio  v  and  Kelvin- Voigt  damping  coefficient  Cd-  As  depicted 
in  Figure  1,  the  axial  direction  is  taken  along  the  x-axis.  The  displacements  of  the  middle 
surface  in  the  longitudinal,  circumferential  and  transverse  directions  are  denoted  by  u,  v  and 
w,  respectively. 

Bonded  to  the  shell  are  s  pairs  of  piezoceramic  patches  which  can  be  employed  as  sensors 
and/or  actuators  in  structural  and  structural  acoustic  applications  [5,  19].  To  simplify  the 
exposition,  the  patches  are  all  asstimed  to  have  thickness  hp^,  Young’s  modulus  jBpe,  Poisson 
ratio  Upe,  and  Kelvin- Voigt  damping  coefficient  cupg.  Furthermore,  it  is  assumed  that  the 
glue  bonding  layer  provides  negligible  contribution  to  the  structural  dynamics.  The  reader  is 
referred  to  to  [6,  9]  for  details  concerning  the  incorporation  of  differing  patch  characteristics 
and  bonding  layers  in  the  ensuing  models). 

Throughout  this  discussion,  it  will  be  assumed  that  external  inputs  to  the  shell  will  be  in 
the  form  of  transverse,  longitudinal  and  circumferential  surface  forces  as  well  cis  line  moments 
and  forces  generated  by  the  patches.  The  surface  forces  can  be  used  to  model  a  variety  of 
phenomena  including  coupling  interactions  with  adjacent  fields  (e.g.,  acoustic  fields  -  see  [5]) 
and  input  from  certain  actuators.  The  patch  moments  and  forces  arise  when  the  elements  are 
used  as  actuators.  While  more  general  inputs  can  be  considered,  the  above-mentioned  cases 
demonstrate  the  flexibility  of  the  numerical  method  for  typical  smart  materials  applications. 


2.1  Strong  Form  of  Equations 

As  detailed  in  [6,  26],  moment  and  force  balancing  yields  the  Donnell- Mushtari  equations 
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as  a  model  for  the  thin  shell  dynamics.  Here  Mx,M0,Mx$,  and  Mox  are  internal  moments, 
while  Nx,  Ne,  Nxo  and  Ngx  denote  internal  force  resultants.  External  surface  forces  are  denoted 
by  qx,qe,Qn  whereas  the  external  resultants  (line  moments  and  forces)  generated  by  the 
patch  pair  are  designated  by  (Mx)pei,  {Me)pei,{Nx)pei  and  {Ne)pei-  The  indicator  function 
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indicates  the  sense  of  the  forces  generated  by  the  patch  pair  having  edges  at  xij,  X2i,  On,  02i. 
The  symmetry  of  the  function  arises  from  the  property  that  for  homogeneous  patches  hav¬ 
ing  uniform  thickness,  equal  but  opposite  strains  are  generated  about  the  point  (xi,  6^  = 
{{xii  + X2i)/2,{0ii  +  02i)/2)  in  the  two  coordinate  directions.  Similarly,  the  characteristic 
function 

.  ..  f  1  ,  Xu<X<  X2i  ,  On  <0  <02i 

=  otherwise 


will  be  used  in  ensuing  discussion  to  isolate  internal  and  external  contributions  due  to  the 
patch  pair. 

Under  the  assumption  that  stress  is  proportional  to  a  linear  combination  of  strain  and 
strain  rate,  the  internal  moments  and  force  resultants  in  regions  of  the  shell  not  covered  by 
patches  are  given  by 
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(see  [6]  or  [26]  for  the  undamped  case).  The  substitution  of  these  equations  in  (2.1)  yields 
the  damped  form  of  the  Donnell-Mushtari  equations  for  a  uniform,  homogeneous  shell  that  is 
devoid  of  patches. 

The  bonding  of  piezoceramic  patches  to  the  shell  produces  contributions  due  to  both 
internal  (material)  and  external  moments  and  forces.  The  internal  contributions  are  due  to 
the  geometrical  and  material  changes  afforded  by  the  patches.  The  external  contributions 
result  from  the  piezoelectric  effect  in  the  patches  which  is  manifested  as  generated  strains  in 
response  to  applied  voltages. 

As  detailed  in  [6,  9],  the  internal  force  resultant  Nx  is  given  by 
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for  systems  in  which  the  patches  have  identical  properties  (compare  with  (2.2)).  Analogous 
expressions  are  given  in  [6,  9]  for  the  remaining  resultants  and  the  resultants  when  the  patches 
have  differing  material  characteristics. 

To  characterize  the  external  contributions,  it  is  typical  to  start  with  the  assumption  that 
the  strains  generated  by  a  patch  are  proportional  to  the  applied  voltage  [6].  Since  differing 
voltages  can  be  applied  to  the  outer  and  inner  patches  in  the  pair,  we  will  differentiate  between 
the  two  with  fji(t)  and  Vi2{t)  used  to  denote  the  voltages  to  the  outer  and  inner  patches  in 
the  pair,  respectively.  The  proportionality  constant  relating  generated  strain  to  the  input 
voltage  is  designated  by  dai.  As  detailed  in  [6],  the  total  external  moments  and  forces  generated 
by  the  patches  are 
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When  substituted  into  (2.1),  the  expressions  (2.4)  provide  the  input  from  the  patches  when 
voltages  are  applied. 
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2.2  Boundary  Conditions 

Appropriate  boundary  conditions  are  dictated  by  the  experimental  setup  or  application  under 
consideration.  In  many  setups,  such  as  the  experimental  shell  apparatus  at  NASA  Langley 
Research  Center,  the  shell  is  supported  by  heavy  endcaps.  In  such  cases,  the  fixed-edge 
conditions 

dw 

u  =  V  =  w  =  =  0  ,  X  =  0,1 

ox 

may  adequately  model  the  end  behavior  of  the  shell. 

If  slight  boundary  movement  or  rotations  are  suspected,  “almost  fixed”  boundary  condi¬ 
tions  of  the  type  discussed  in  [7,  8]  can  be  employed.  Such  boundary  conditions  can  be  used 
to  model  the  slight  energy  which  result  due  to  imperfect  clamping  of  the  structure. 

A  third  type  of  commonly  considered  boundary  conditions  are  the  simply-supported  or 
shear  diaphragm  edge  conditions 

V  =  w  =  Mx  =  Nx  =  0  ,  x  =  0,£.  (2-7) 

These  boundary  conditions  are  theoretically  attractive  since  they  provide  one  of  the  few  sit¬ 
uations  in  which  explicit  modal  expansions  can  be  calculated  for  the  Donnell-Mushtari  shell 
equations.  They  are  of  limited  use  when  modeling  experimental  shells,  however,  since  they 
are  appropriate  only  for  endcaps  which  prevent  movement  in  the  v  and  w  directions  but 
are  sufficiently  flexible  in  the  x-direction  so  that  negligible  moments  Mx  and  forces  Nx  are 
generated. 

The  Galerkin  method  of  this  work  is  equally  efficient  to  implement  for  models  incorporating 
the  fixed-edge  boundary  conditions  (2.6),  “almost  fixed”  boundary  conditions,  or  the  simple- 
supported  edge  conditions  (2.7).  The  method  is  demonstrated  for  the  fixed-edge  and  simply- 
supported  conditions  while  modifications  to  adapt  the  method  to  “almost  fixed”  can  be  found 
in  [7,  8].  Later  discussion  will  also  illustrate  the  manner  through  which  the  method  can  be 
adapted  to  alternative  boundary  conditions  which  many  arise  in  physical  applications. 

2.3  Weak  Form  of  Equations 

As  noted  in  (2.1),  the  strong  form  of  the  modeling  equations  involves  first  and  second  deriva¬ 
tives  of  both  the  moment  and  force  resultants.  For  structures  with  surface-mounted  or  embed¬ 
ded  actuators  or  sensors  (e.g.,  piezoceramic  patches),  this  differentiation  leads  to  difficulties 
due  to  discontinuities  in  the  resultants.  The  internal  resultants  contain  piecewise  discontinu¬ 
ities  due  to  material  changes  introduced  by  the  patches  (see  for  example,  (2.3)).  The  external 
patch  contributions  are  also  discontinuous  since  they  are  defined  only  over  regions  covered  by 
patches  (see  (2.4)).  As  a  result,  formal  analysis  and  approximation  using  the  strong  form  of 
the  modeling  equations  leads  to  derivatives  of  the  Dirac  delta  “function.” 

To  alleviate  these  difficulties,  we  consider  a  weak  form  of  the  modeling  equations.  As 
detailed  in  [6],  such  equations  can  be  derived  from  energy  considerations  and  are  equivalent  to 
the  strong  form  of  the  modeling  equations  under  suitable  smoothness  conditions.  This  yields 
a  form  of  the  model  which  facilitates  well-posedness  analysis  [5]  and  eliminates  the  difficulties 
associated  with  the  discontinuous  resultants.  Moreover,  smoothness  requirements  on  basis 
functions  are  reduced  which  proves  advantageous  when  constructing  system  matrices. 


7 


The  weak  form  of  (2.1),  as  derived  from  energy  considerations  in  [6],  is  given  by 

r  c  ° 


(2.8) 


|^2■!T  ri  I  0‘-yj 

I  l\Rph^m  +  N.r„-RM. 


d^w 


dx^  R  "  de^ 


■m.’is + E  \m^u^ + 

2=1  L 


) 


dxdd 


dxdO  =  0 


for  all  rf  —  [rj\^rj2,rj-^  €  V  where  V  denotes  the  space  of  test  functions  (as  detailed  in  [6],  the 
indicator  functions  appear  only  in  the  definition  of  the  external  force  resultants  in  the  weak 
form  of  the  equations).  Specific  choices  for  V  depend  upon  the  boundary  conditions.  For  the 
fixed-edge  boundary  conditions  (2.6),  the  space  of  test  functions  is  taken  to  be 

V  =  H^{n)  X  H^{n)  X  H^(n) 


where 

Hlitt)  =  {,  e  /?'(«)  :  >)(0)  =  m  =  0} 

HUa)  =  {<,e  H^CI) :  7,(0)  =  ,.(o)  =  =  ri4e)  =  0}  . 

In  the  case  of  the  simply-supported  edge  conditions  (2.7),  essential  boundary  conditions  are 
imposed  only  on  the  circumferential  and  transverse  functions  and  V  is  taken  to  be 

F  =  H^{n)  X  Hl{n)  X  Hl(n) 


with 

Hliil)  =  {7,  6  H\n)  :  7,(0)  =  r,(e)  =  o}  . 

In  general,  V  is  simply  taken  as  the  subset  of  the  traditional  Sobolev  spaces  satisfying  essential 
boundary  conditions  (the  reader  is  referred  to  [7]  for  modifications  to  account  for  “almost 
fixed”  boundary  conditions). 

A  comparison  between  (2.8)  and  (2.1)  illustrates  that  in  the  weak  form,  derivatives  are 
transferred  from  the  discontinuous  resultants  onto  suitably  smooth  test  functions.  This  allevi¬ 
ates  the  difficulties  associated  with  the  discontinuities  and  reduces  smoothness  requirements 
on  approximate  solutions. 
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3  Modal  Solution  —  Simply- Supported  Boundary  Con¬ 
ditions,  Constant  Coefficients,  No  Damping 


As  noted  in  the  introduction,  closed  form  modal  expressions  for  shell  models  can  be  obtained 
only  in  a  limited  number  of  ceises.  One  case  in  which  analytic  frequencies  and  modes  can 
be  calculated  is  the  undamped  (c^  =  0)  shell  model  with  simply-supported  edge  conditions 
(2.7)  and  constant  parameters  p,  v  and  E.  Separation  of  variables  in  this  setting  is  classical 
and  can  be  found  in  numerous  texts  (e.g.,  [27]);  we  summarize  the  arguments  here  to  facili¬ 
tate  numerical  examples  in  later  sections.  While  such  conditions  are  not  attained  in  typical 
shell  applications,  the  consideration  of  the  undamped,  constant  coefficient  shell  model  with 
simply-supported  edge  conditions  provides  an  excellent  means  of  testing  discretization  tech¬ 
niques  since  approximate  solutions  can  be  compared  with  analytic  values.  The  discretization 
techniques  can  then  be  used  to  approximate  natural  frequencies,  modes  and  shell  dynamics  in 
models  which  incorporate  piecewise  constant  parameters  (including  damping)  and  physically 
realistic  boundary  conditions. 

Throughout  this  section,  it  is  assumed  that  c^  =  0  and  p,  v  and  E  are  constant.  Since 
our  interest  here  is  restricted  to  the  calculation  of  natural  frequencies  and  modes,  we  will  also 
consider  the  shell  model  to  be  unforced  (no  external  patch  contributions  and  qx  =  qe  =  =  ^)- 

In  this  case,  the  strong  form  of  the  Donnell-Mushtari  thin  shell  equations  (2.1)  can  be  written 
in  the  operator  format 


=  Lu 


(3.1) 


where  u  =  [u,  u,  wY  ■  The  operator  L  here  is  given  by 


Eh 


(1  -  v) 

dx^  ‘IE?  de^ 
(1-1- 1/) 


2R  dxdO 


2R  dxde 

(l-i/)  92  1  92 

2  dx^  R?  902 


V  9 
Rdx 

J__9 

R?de 


1/  9 

1 

9 

_ _ ji 

94 

2 

94 

1 

94  ■ 

'Wdx 

90 

R?  * 

9x4  ^ 

R? 

^  ie4 

904 

where  k  =  h?l  12.  Boundary  conditions  for  the  shell  are  denoted  by 

Bu  =  0  . 

Due  to  linearity,  the  displacements  are  expressed  as 

u(t,  x,  9)  =  T{t)U{x,  6) 
v(t,x,d)  =  T{t)V{x,e) 
w(t,x,9)  =  T{t)W(x,9)  , 


and  spatial  and  temporal  components  are  separated  to  yield  the  eigenvalue  problem 


along  with  the  temporal  equation  T"  +  =  0.  Here  u  is  the  circular  frequency  of  vibration 

and  U  =  [U,  V,  W]^  contains  the  vibration  modes  in  the  axial,  circumferential  and  transverse 
directions.  We  point  out  that  separating  variables  to  arrive  at  the  eigenvalue  problem  (3.2)  is 
equivalent  to  the  assumption  of  a  harmonic  response  in  all  components  and  one  obtains  the 
same  eigenvalue  problem  in  both  cases. 

The  structure  of  the  eigenfunctions  U,  V  and  W  is  dependent  upon  the  boundary  conditions 
with  closed  form  expressions  attainable  in  only  a  few  cases.  To  illustrate,  we  consider  the 
eigenvalue  problem  (3.2)  with  the  simply- supported  boundary  conditions  (2.7).  One  form  of 
the  natural  vibration  modes  for  this  case  is 


(3.3) 


/  rr'TTX  \ 

U{x,0)  =  Ai  cos  ( ]  cos(m^) 
V{x,6)  —  Bi  sin  sm{m6) 

(TITTX  \ 

— cos{m0)  . 


The  choice  of  the  cosine  representation  for  the  x  component  of  the  axial  vibration  and  sines 
for  the  circumferential  and  radial  vibrations  guarantees  that  the  modes  satisfy  the  boundary 
conditions.  The  relationship  between  the  circumferential  and  axial/radial  vibrations  dictates 
that  the  former  must  be  90°  out-of-phase  from  the  latter  two  in  6  (see  the  moment  and  force 
expressions  (2.2)  or  operator  definition  (3.1)). 

To  determine  the  relationship  between  the  frequencies  oj  and  the  wave  numbers  m  and  n, 
the  expressions  (3.3)  are  substituted  into  (3.2)  to  yield  the  system 


where 


■  uX  ' 

■  ■ 

■  0  ' 

^(1  Av)  -  H2  -m 

== 

0 

uX  —m  —  H3 

.  C'l  . 

.0  _ 

H2  =  y(1  -  z/)  -f 
H3  =  l  +  k(^X^  + 


(3.4) 


and 


P_ 
E 
nirR 


=  .^(1  _ 
E 

x  = 


I-  ^ 

~  R?~  12R?  ' 

A  nontrivial  solution  is  determined  by  setting  the  determinant  to  zero.  This  yields  the 
cubic  equation 

-f  KrO?  -  ATo  =  0  (3.5) 
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in  0,^.  The  coefficients  here  are  given  by 
K2  =  1  +  ^(3  -  u) 

Ki  =  ^(1  -  v)  ^(3  +  2i/)A^  +  m^  +  (m^  +  A^)  +  j 

^0  =  ^(1  -  i')  [(1  -  +  k  (m^  +  A^)^j  . 

The  solutions  to  this  cubic  equation  are  then  given  by 

Km  =  -yKl  -  3^1  cos  ~  j  ^  ,  j  =  1.2,3 


where 


a  =  cos 


-1 


-21  Kq  -  2Kl  +  ^K^K2 
2^J{Kl  -ZK^f 


Finally,  the  natural  frequencies  for  the  shell,  with  units  of  hertz,  are  given  by 


r  _  flmnj 

Jmnj  — 


E 


i  =  1,2,3  . 


(3.6) 


27rR\l  p{l-u^) 

We  point  out  that  due  to  the  cubic  nature  of  the  characteristic  equation,  three  natural  fre¬ 
quencies  and  mode  shapes  are  obtained  for  each  set  of  wave  numbers  m  and  n.  This  leads  to 
a  significant  interlacing  of  frequencies. 

A  second  set  of  modes  which  satisfy  the  boundary  conditions  and  6  compatibility  criteria 
is 

U{x,  6)  =  A2  cos  sin(m0) 


V (x,  0)  =  B2  sin  f  )  cos{m9) 

W{x,6)  =  C^jsin  sm(m0) 


(3.7) 


The  system  which  arises  when  using  this  form  of  the  modes  is 


1 

G 

1 _ 

-^(1  +  ^)  ^A  ■ 

A2 

■  0  ■ 

n^-ff2  m 

B2 

= 

0 

I/A 

m 

.  C2  . 

0 

which  differs  in  sign  in  the  (1, 2),  (2, 1),  (2, 3)  and  (3, 2)  elements  from  (3.4)  which  was  obtained 
using  the  first  set  of  modal  expansions.  The  characteristic  equations,  however,  are  the  same 
in  both  cases  and  are  given  by  (3.5).  Hence  both  sets  of  expressions  yield  the  same  natural 
frequencies  for  the  shell.  This  should  be  expected  since  the  second  set  simply  represents  a 
phase  shift  in  0  of  the  first.  For  this  reason,  most  authors  consider  only  the  first  set  of  vibration 
modes. 

In  some  instances,  however,  the  second  set  (3.7)  contributes  linearly  independent  modes 
and  hence  must  be  retained  in  order  to  obtain  a  complete  bcisis  for  approximation.  This  is 
illustrated  in  the  examples  of  Section  5.1.  Further  details  regarding  the  modal  analysis  can 
be  found  in  [16]. 
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4  Approximation  Method 

To  approximate  the  solutions  u,v,w  to  the  system  (2.8)  or  eigenvalue  problem  (3.2),  basis 
functions  and  B.uj^{0,x)  in  V  are  chosen  and  used  to  form  expansions 

u^(t,  x)  =  ^  Uk{t)Buk{9,  x) 

k=i 

v^{t,9,x)  =  '^Vk{t)B^^{e,x)  (4.1) 

k=l 

Afw 

w^{t,9,x)  =  Y^Wkit)By;^{9,x)  . 

k=i 


The  time-dependent  generalized  Fourier  coefficients  Uk{t)^Vk{t),Wk{t)  are  then  determined  by 
orthogonalizing  the  residual  with  respect  to  linearly  independent  test  functions  from  the  span 
of  the  bases. 

To  exploit  the  tensor  nature  of  the  cylindrical  shell  domain  fl  and  periodicity  in  9,  the 
bases  for  the  three  displacements  are  constructed  with  Fourier  components  in  9  and  spline 
components  in  x.  When  expressed  with  real  components,  this  yields 


u^{t,  9,  x) 


+ 


Nu 

Y^UOn{t)Bu„{x) 

n=l 

Mu  Nu  Mu  Nu 

m=l  n=l  771=1 71=1 


(4.2) 


with  similar  expressions  for  v^{t,  9,  x)  and  w^{t^  9,  x).  In  the  expansions,  M^,  are  the 

Fourier  limits  and  Nu,  Nu,  denote  the  number  of  splines  Bun{x),  jB„„(x),  B^uix)  used  when 
approximating  the  longitudinal,  circumferential  and  transverse  displacements,  respectively. 
The  total  number  of  basis  functions  for  each  component  is  Nu  =  Nu  •  (2Mu  +  l),W„  =  Ny  • 
(2Mu  +  l)  and  AC  =  Nw  ■  (2AC  +  !)• 


4.1  Axial  Basis  Functions 

The  choice  of  splines  in  the  axial  variable  is  motivated  by  the  following  criteria: 

(i)  efficiency; 

(ii)  flexibility  with  regard  to  internal  and  external  patch  contributions; 

(iii)  adaptability  with  respect  to  various  boundary  conditions; 

(iv)  accuracy. 

For  the  longitudinal  and  circumferential  displacements,  both  linear  and  cubic  spline  bases  are 
considered.  Due  to  differentiability  requirements,  only  cubic  splines  are  used  when  approxi¬ 
mating  the  transverse  displacement. 

In  all  cases,  a  uniform  partition  along  the  x-axis  is  considered  with  gridpoints  x„  = 
nhx  ,  hx  —  IJN  ,  n  =  1,  •  •  • ,  A^.  For  n  =  — 1,0, 1,---,W-|-1,  standard  cubic  splines  are  de¬ 
fined  by 
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^n(^)  — 


hi 


f  2) 

"i"  Sh^(^x  Xjj— 1)  “I"  3hx{^x  Xfi—i)  8(2^  Xn—2') 

hx  “I"  2;)  -|-  3/ix(^n+i  3(2;ji^.i  x) 


(a:n+2 

0 


0^ 


X  e  [x„_2,X„_i] 

X  €  [Xn_i,X„] 

a:  €  [Xn,  2:^+1] 

X  G  [a^n+i)  a;^^2] 
otherwise 


(see  [29]).  The  standard  linear  splines  are  defined  by 

(  (x  Xji—l)  )  2;  £  [Xji_i,Xn] 


Cn{x)  = 


h. 


(Xn+i  -  X)  ,  X  G  [x 

n?  a^n+i] 

0  ,  otherwise 


(4.3) 


(4.4) 


bo{x)  —  46_i(x) 

,  n  =  0 

6i(x)  -  6_i(x) 

,  n  =  1 

bn{x^ 

,  n  =  2,  •••, 

bN-\{x)  —  bN+l 

,  n  =  N  —  1 

bN{x)  —  46jv+i(a:) 

,  n  =  N 

for  n  =  0,  •  •  •  ,iV. 

In  both  cases,  the  standard  splines  must  be  modified  to  satisfy  the  essential  boundary 
conditions  to  ensure  approximate  solutions  in  V.  These  modifications  are  summarized  as 
follows. 

(i)  Fixed  displacements  at  x  =  0,^: 

To  construct  basis  functions  satisfying  fixed  displacements  (but  unspecified  slopes)  at 
X  =  0,  the  modified  cubic  splines  are  taken  to  be 


(4.6) 


for  a  total  of  iV  +  1  functions.  To  satisfy  the  same  condition,  the  first  and  last  linear  splines 
are  omitted  from  the  set  to  yield 

4(x)  =  Cn(x)  ,  n  =  I,--- ,iV- 1  .  (4.6) 

It  should  be  noted  that  with  these  definitions,  6n(0)  =  hn{l)  =  0  and  c,i(0)  =  Cn(^)  =  0. 

(ii)  Fixed  displacements  and  slope  at  x  =  0,£: 

Only  the  cubic  splines  are  required  to  satisfy  a  condition  of  fixed  displacement  and  slope 
at  the  endpoints  since  this  is  a  condition  imposed  only  on  the  transverse  displacements.  In 
this  case,  the  modified  splines  are  taken  to  be 

bo(x)  —  26_i(x)  —  26i(x)  ,  n  =  1 

bn{x)=<  5„(x)  ,  n  =  2,  •••,iV-2  ,  (4.7) 

_  6^-(x)  -  2&^-_i(x)  -  2bN+i{x)  ,  n  =  N  -1 

for  a  total  of  iV  —  1  basis  functions.  Note  that  these  functions  satisfy 

640)  =  ~6;(0)  =  k{£)  =  I'M  =  0  . 
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4.2  Bases  for  Specific  Boundary  Conditions 

Appropriate  bases  for  various  boundary  conditions  are  then  constructed  by  considering  mod¬ 
ified  splines  satisfying  any  essential  boundary  conditions.  For  example,  the  bases  used  for 
simply-supported  shells  (2.7)  must  satisfy  the  conditions  jBu„(0)  =  =  0,  = 

(^)  ==  0  (the  moment  and  shear  conditions  are  natural  and  hence  do  not  need  to  be  explic¬ 
itly  enforced).  Choices  for  shell  models  with  fixed-edge  conditions  are  summarized  in  Table  1 
while  corresponding  choices  for  simply-supported  edge  conditions  are  given  in  Table  2.  The 
Icist  column  in  each  table  summarizes  the  total  number  of  axial  functions  in  each  expansion. 

Once  the  axial  bases  are  chosen,  they  are  combined  with  the  Fourier  components  in  6 
to  form  the  total  bases.  We  reiterate  that  the  basis  limits  in  each  case  are  given  by  A4  = 
Nu  '  (2Afu  -|- 1)  ,A/^  =  Ny  •  {2My  -b  1)  and  Afy,  =  Ny,  •  {2Myj  -|- 1). 

We  point  out  that  the  linear  and  cubic  spline  bases  described  here  are  but  two  choices  from 
among  many  that  can  be  made  for  the  axial  components.  For  the  applications  considered  here, 
cubic  splines  provided  a  good  balance  between  accuracy,  eflSciency  and  adaptability  with  re¬ 
gards  to  patches  and  boundary  conditions.  If  higher  accuracy  is  desired,  however  (with  slightly 
more  expense  when  constructing  system  matrices),  quintic  splines  can  be  employed.  Similarly, 
spectral  expansions  employing  Legendre,  Chebyshev  or  sine  functions  can  be  employed  once 
modifications  have  been  made  for  boundary  conditions. 


Shell  Component 

Axial  Basis  Functions 

Component  Definition 

Axial  Limit 

longitudinal  -  linear 

(^)  ~ 

Cn{x)  defined  in  (4.6) 

II 

1 

cubic 

Bu„{x)  =  bn{x) 

bn{x)  defined  in  (4.5) 

Ny,  =  Ny  +  1 

circumferential  -  linear 

Bv„{x)  =  Cn{x) 

Cn(x)  defined  in  (4.6) 

II 

1 

cubic 

<rJ^ 

II 

C 

bn(x)  defined  in  (4.5) 

Ny  -  Ny  +  1 

transverse  -  cubic 

—  ^n(^) 

bn(x)  defined  in  (4.7) 

1 

II 

Table  1.  Axial  basis  definitions  for  models  with  fixed-edge  boundary  conditions  (2.6). 


Shell  Component 

Axial  Basis  Functions 

Component  Definition 

Axial  Limit 

longitudinal  -  linear 

Buni^)  = 

Cn{x)  defined  in  (4.4) 

Ny  =  Ny  +  1 

cubic 

bn{x)  defined  in  (4.3) 

Ny  =  Ny  +  3 

circumferential  -  linear 

Cn{x)  defined  in  (4.6) 

Ny=Ny-l 

cubic 

A 

bn{x)  defined  in  (4.5) 

Ny  =  Ny  -f-  1 

transverse  -  cubic 

—  ^71  (^) 

A 

bn{x)  defined  in  (4.5) 

II 

-f 

Table  2.  Axial  basis  definitions  for  models  with  simply-supported  boundary  conditions  (2.7). 
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4.3  Components  in  the  Finite  Dimensional  System 

With  bases  thus  defined,  approximating  subspaces  are  taken  to  be  =  span{B„  = 

span{;Bi,^};t=i  h  should  be  noted  that  x  x  C  V. 

The  approximating  system  is  then  determined  by  restricting  the  weak  form  (2.8)  to  with 
basis  functions  used  as  test  functions.  This  is  equivalent  to  orthogonalizing  the  residual  with 
respect  to  elements  from  . 


4.4  Matrix  System 

To  form  the  resulting  matrix  system,  the  generalized  Fourier  coefficients  for  the  three  expan¬ 
sions  (4.1)  are  first  consolidated  in  the  vectors 


Ui{t) 

,  V^^(t)  = 

Ul(t) 

,  w-^“'(t)  = 

Wi{t) 

.  WAr„(t)  . 

.  W.(t)  . 

(4.8) 


The  full  set  of  coefficients  is  then  represented  as  •d^(i)  =  [^(t)''^“,  V"^‘'(t),  W'^“'(t)]^,  where 

.A/”  =  J\fu  +  M-u  d-  A/^. 

The  mass,  stiffness  and  damping  matrices  as  well  as  the  forcing  vector  for  the  full  system 
have  the  form 


'  Um 

= 

Vm 

Wm  _ 

and 


Fix  +  Fi2 

kFii 

Ki  = 

U2I  +  U22 

V21  +  ^22 

W21 

U31 

F31 

- 1 

tH 

II 

W 

’  ^11  +  Ui2 

+  Fi2 

1 - 

rH 

II 

U2I  +  U22 

+  V22 

W21 

Usi 

V31 

ELi  M's).  . 

(4.9) 


(4.10) 


The  various  submatrices  contain  individual  components  which  arise  when  the  weak  form  (2.8) 
is  restricted  to  .  For  example,  the  approximation  of  the  mass  and  stiffness  components  in 
the  longitudinal  equation  of  (2.8)  yields 
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r2'iT  r£  I  ^ 

(i)  WM]j,k  =  ^ph+2j2PpehpeXpei  RB^^BujdxdO 

dBu,  dBu, 


nt 


EhR  2EpehpeR 

Y~;^  +  ^  1  _  1,2 


i=l 


pe 


dx  dx 


dxdO 


nt 

_ 


Ehi/ 


+  ^- 


^Eryehrye^r) 


nt 

. 


1  —  i/2  ‘  ^  \  _  yl 

t—\  ^  pe 


Ehu  — \  2E‘Qeh 


'Xpe,- 


9a; 


dxdO 


\^y2 


+  ^-T 


Z=1 


'pe**'pe^pe 
pe 


dBu 

B^,-^dxde 


r  Eh  ^  f 

WT7)^?ri 


Rpe^pe 


r27T  r£ 

(vi)  |^2iJ(l  +  1.)  '  ^  il(l  +  Z/pe) 

/•2ir  f  s  "I 

(vii)  [Fu]j  =  l  J^^Rq.Bu,+Rj2i^-)pei-^jdxd0 

with,  similar  expressions  for  the  submatrices  f/n,  l^ii?  1^2  and  Uu  which  contain  the  in¬ 
ternal  damping  contributions.  The  submatrices  for  the  approximated  circumferential  and 
transverse  equation  are  constructed  analogously  with  details  given  in  [16].  It  should  be  noted 
that  due  to  the  presence  of  the  characteristic  functions  in  the  force  and  moment  resultants, 
the  domains  of  integration  for  the  patch  contributions  are  restricted  to  those  regions  covered 
by  patches.  This  is  consistent  with  the  fact  that  material  and  geometric  changes,  as  well  as 
patch  inputs,  occur  only  in  those  regions. 

In  first-order  form,  the  system  has  the  form 


2(1  +  *')  ^  1  +  V 

Eh 


Xpei 


dxdO 


Epehpe 


dBy,  dBu, 

dx  do 

dBu,  OBuj 


7Xpei 


do  do 


dxdO 


'Kg 

0 

■  d^{t)  ■ 

0 

.  . 

_ 

' 

0 

■  t?^(o)  ■ 

r 

0 

_ 

_  i?^(0)  _ 

0 

■Kg 


Kg 

-Kg 


>>^(*)  J 


+ 


0 

F^{t) 


91 

92 


(4.11) 


Multiplication  of  the  inverted  mass  matrix  yields  a  Cauchy  equation  of  the  form 

y^{t)  =  AV{t)+9^{t) 

2/^(0)  =  Vo  ) 


(4.12) 
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where  G  .  This  yields  a  form  which  is  suitable  for  simulations, 

parameter  estimation  and  control  applications.  Note  that  the  system  can  also  be  adapted  to 
alternative  boundary  conditions  through  modifications  of  the  first  and  last  basis  functions. 
Flexibility  in  this  regard  is  also  a  hallmark  of  the  method. 

Finally,  we  point  out  that  approximation  of  the  eigenvalue  problem  (3.2)  using  the  Fourier- 
Galerkin  method  yields  the  generalized  matrix  eigenvalue  problem 

(4.13) 

where  and  are  defined  in  (4.9).  As  illustrated  in  the  examples  of  Section  5.1,  (4.13) 
can  be  used  to  approximate  the  frequencies  and  modes  for  undamped  shells  with  various 
internal  characteristics  (due  to  actuators  and  sensors)  and  boundary  conditions. 

5  Examples 

To  illustrate  various  facets  of  the  Fourier-Galerkin  method,  two  types  of  examples  are  pre¬ 
sented.  In  the  first  set  of  examples,  approximation  of  frequencies  and  modes  is  considered 
for  simply-supported  and  fixed-edge  shells.  This  is  accomplished  by  forming  the  mass  and 
stiffness  matrices  and  solving  the  generalized  matrix  eigenvalue  problem  (4.13).  This  serves 
two  purposes.  It  illustrates  the  accuracy  of  the  method  in  approximating  the  mass  and  stiff¬ 
ness  components  in  the  system  and  demonstrates  the  method  as  an  effective  technique  for 
obtaining  modal  information  when  analytic  or  experimental  values  are  unavailable. 

The  second  set  of  examples  illustrates  the  approximation  of  static  and  dynamic  shell 
responses  to  various  forces.  The  initial  steady  state  example  illustrates  that  the  expected 
0{h^)  and  0{h^)  convergence  rates  are  obtained  when  approximating  a  static  shell  response 
(this  example  involves  the  solution  of  a  matrix  system  with  coefficient  matrix  and  vector 
F-^).  Analogous  dynamic  examples  demonstrate  that  the  same  accuracy  is  obtained  when  the 
method  is  used  to  obtain  the  time-dependent  system  (4.12)  which  is  then  marched  in  time. 
The  latter  example  also  incorporates  the  internal  Kelvin- Voigt  damping.  The  dynamics  of  the 
system  in  response  to  patch  inputs  are  considered  in  the  final  example.  Natural  frequencies  are 
compared  to  those  obtained  by  solving  the  generalized  matrix  eigenvalue  problem  to  illustrate 
consistency  among  techniques. 

In  these  examples,  the  following  shell  characteristics  were  used:  thickness  h  =  .01  in, 
length  ^  =  12  in,  radius  R  =  3.0  in,  mass  density  p  =  .283  Ibfin^,  Young’s  modulus  =  3.0  x 
10^  Ibfin^,  Poisson  ratio  p  =  0.3  and  Kelvin- Voigt  damping  coefficient  cd  =  15.09936  Ib-in-s^. 
These  values  were  chosen  so  as  to  permit  comparison  with  Donnell-Mushtari  shell  results  in 
[24,  pp.  307-310].  We  emphasize  that  the  method  is  not  restricted  to  such  dimensions  or 
ratios,  however,  and  analogous  convergence  results  have  been  obtained  with  a  wide  range  of 
parameters. 

5.1  Modal  Examples:  Simply- Supported  Shell 

In  Section  3,  separation  of  variables  was  used  to  obtain  analytic  expressions  for  modes  and 
natural  frequencies  of  shell  models  with  simply-supported  boundary  conditions.  Here  we  com¬ 
pare  the  approximate  solutions  obtained  by  solving  the  matrix  eigenvalue  problem  (4.13)  with 
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the  analytic  values.  This  provides  a  means  of  testing  the  accuracy  of  the  method  before  using 
it  to  approximate  in  settings  in  which  analytic  solutions  are  unavailable.  To  illustrate,  we 
consider  the  cases  of  purely  axisymmetric  modes  (m  =  0,  n  >  1),  purely  extensional  modes 
(m  >  l,n  =  0)  and  general  modes  (m  >  l,n  >  1). 

5.1.1  m  =  0,  n  >  1  Axisymmetric  Modes 

Analytic  frequencies  for  the  axisymmetric  case  can  be  calculated  via  (3.6)  with  m  =  0  while 
corresponding  modes  can  be  determined  from  (3.3)  or  (3.7)  (see  [16]  for  details).  Corre¬ 
sponding  approximate  frequencies,  obtained  by  solving  (4.13)  with  Fourier  limit  M  =  0,  are 
compared  with  analytic  values  in  Table  3.  As  noted  in  Section  3,  three  frequencies  are  ob¬ 
tained  for  each  configuration  of  1/2  wavelength  (fixed  values  of  n),  and  this  is  reflected  in  the 
analytic  and  approximate  results  in  Table  3.  Two  of  the  frequencies  correspond  to  coupled 
axial/radial  modes  while  the  remaining  one  corresponds  to  an  uncoupled  torsional  (circum¬ 
ferential)  mode.  The  frequencies  of  the  torsional  modes  are  boxed  to  facilitate  comparison 
between  the  corresponding  functional  modal  approximations. 

The  cubic  spline  approximates  in  Table  3  were  obtained  using  N  =  S  basis  functions 
whereas  iV  =  32  linear  splines  were  used  to  obtain  the  corresponding  results  in  columns  8-10. 
It  is  noted  that  accurate  frequency  approximates  are  obtained  with  the  cubic  splines  with 
a  maximum  relative  error  (over  the  reported  results)  of  0.3%  occurring  for  the  frequency 
2255.52  Hz.  Due  to  the  limited  accuracy  of  the  linear  splines,  significantly  more  basis  functions 
must  be  used  to  obtain  accurate  approximations,  and  a  1.25%  error  at  2255.52  Hz  remains, 
even  with  N  =  32  functions. 

For  n  =  1,  cross  sections  at  0  =  0  of  corresponding  modes  approximated  using  M  =  0,N  = 
16  cubic  splines  are  plotted  in  Figure  2.  The  coupling  between  the  axial/radial  components 
is  readily  verified  for  the  406.8  Hz  and  603.83  Hz  modes  while  the  plot  of  the  266.05  Hz 
mode  illustrates  that  for  M  =  0,  the  torsional  mode  is  uncoupled  (examples  in  Section  5.1.3 
illustrate  that  all  three  modes  are  coupled  for  M  >1). 

It  should  be  noted  that  for  each  n,  all  three  frequencies  and  modes  are  automatically 
yielded  by  the  approximation  method  whereas  (3.3)  and  (3.7)  must  be  employed  to  obtain 
analytic  expressions  for  the  axial/radial  and  torsional  modes,  respectively.  In  this  aspect,  the 
approximation  method  facilitates  the  calculation  of  a  complete  modal  set. 


n 

Analytic  Frequencies 

Cubic  Spline  Approx 

1 

266.05 

406.80 

406.80 

603.83 

2 

531.53 

531.53 

532.11 

924.29 

3 

Kami 

1362.11 

798.22 

1362.19 

4 

543.53 

1064.22 

1807.86 

543.58 

1064.82 

1808.81 

5 

544.60 

1330.27 

2255.52 

544.85 

1333.95 

2262.44 

Linear  Spline  Approx 


603.87 

925.95 

1368.07 

1822.21 

2283.78 


Table  3.  Analytic  frequencies  and  approximate  values  obtained  using  N  =  8  cubic  splines 
and  N  —  32  linear  splines.  Frequencies  of  the  torsional  modes  are  boxed. 
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U^{x,Q)  F^(x,0)  W^{x,Q) 


1330.27  Hz 

C^^(x,0)  F^(x,0)  T7^(x,0) 


2255.52  Hz 

Figure  2.  Approximate axisymmetric modes  U^{x,Q),V'^{x^Qi)  and  W''^(x,0)  corresponding 
to  the  544.60,  1330.27  and  2255.52  Hz  natural  frequencies. 
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5.1.2  m  >  l,n  =  0  Purely  Extensional  Modes 

A  second  case  in  which  analytic  expressions  for  the  shell  frequencies  and  modes  are  easily  ex¬ 
pressed  occurs  when  purely  extensional  modes  are  present.  The  first  three  analytic  frequencies 
given  by  (3.6)  are  compared  in  Table  4  with  approximate  values  obtained  with  N  =  S  cubic 
splines.  Due  to  the  accuracy  of  the  method,  the  approximate  frequencies  agree  with  analytic 
values  to  within  at  least  two  decimal  places  in  all  three  cases.  Corresponding  modal  plots  are 
given  in  Figure  3.  The  numerical  plots  in  the  left  column  illustrate  the  two  approximate  solu¬ 
tions,  for  each  frequency,  obtained  by  solving  the  matrix  eigenvalue  problem  (4.13)  (again,  the 
full  solution  set  is  automatically  obtained  by  the  approximation  method).  The  corresponding 
analytic  modes  given  by  (3.3)  and  (3.7)  are  depicted  in  the  right  column  of  the  figure.  With 
n  =  0,  it  is  seen  from  (3.3)  and  (3.7)  that  the  analytic  expressions  for  the  longitudinal  modes 
are  U(x,9)  =  Ai  cos(m0)  and  U{x,9)  —  A2s\n{m9)  with  zero  displacement  expected  in  the 
circumferential  and  transverse  components.  This  longitudinal  behavior  is  noted  for  both  the 
analytic  and  approximate  modes  in  Figure  3. 


Analtyic  (Hz)  Galerkin  (Hz) 

m  =  l,n  =  0 
m  =  2,n  =  0 
m  =  3,  n  =  0 

338.75  338.75 

677.50  677.50 

1016.25  1016.25 

Table  4.  Analytic  and  approximate  frequencies  for  purely  extensional  modes  obtained  using 
N  =  S  cubic  splines. 

5.1.3  m  >  l,n  >  1  General  Shell  Modes 

For  the  general  case  with  m  >  1,  n  >  1,  the  axial,  radial  and  torsional  modes  are  all  coupled 
with  analytic  values  for  the  natural  frequencies  determined  by  (3.6)  where  one  of  the 

three  solutions  to  the  cubic  equation  (3.5).  These  analytic  values  are  compared  in  Table  5 
with  approximate  frequencies  obtained  using  the  cubic  spline  basis  with  N  =  S,M  =  3.  For 
comparison  sake,  frequencies  of  the  axisymmetric  (m  =  0)  and  purely  extensional  (n  =  0)  are 
also  included  in  this  table.  The  accuracy  of  the  cubic  spline  discretization  leads  to  approximate 
frequencies  having  relative  errors  less  than  0.5%  for  the  reported  values  (the  largest  relative 
error  for  the  reported  values  occurs  in  the  second  frequency  obtained  with  m  —  3,  n  =  5).  The 
convergence  of  the  method  and  accuracy  obtained  using  the  cubic  splines  is  further  illustrated 
by  the  frequency  results  in  Table  6  where  approximates  obtained  with  N  —  8  and  iV  =  16  are 
compared.  A  check  of  the  relative  errors  shows  that  the  method  is  converging  more  quickly 
than  the  expected  0{h^)  rate. 

The  modes  for  this  general  case  are  again  obtained  from  (3.3)  and  (3.7).  The  coupling 
between  the  axial,  radial  and  torsional  modes  can  be  seen  in  Figure  4  where  the  873.61  Hz 
mode  is  plotted  along  the  axial  line  £  —  {(x,  9)\9  =  7r/4, 0  <  x  <  12}. 

We  reiterate  that  the  full  set  of  frequencies  and  modes  can  be  approximated  using  the 
Fourier- Galerkin  method  by  simply  solving  the  matrix  eigenvalue  problem  (4.13).  Moreover, 
this  technique  can  be  directly  extended  to  the  problem  of  calculating  frequencies  and  modes  for 
shells  with  other  boundary  conditions  by  suitably  modifying  the  spline  basis.  This  is  illustrated 
in  the  next  section  where  the  modal  analysis  of  a  shell  with  clamped  ends  is  considered. 
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338.75  Hz 


677.50  Hz 


1016.25  Hz 


Figure  3.  Purely  extensional  modes  in  tiie  longitudinal  direction.  For  N  =  S  cubic  splines, 
numerical  values  of  along  the  curve  i  —  {{x,  9)\0  <  9  <27: ,  x  —  0}  are  plotted  in  left  column 
Corresponding  shell  modes  given  by  (3.3)  and  (3.7)  with  n  =  0  are  depicted  in  right  column 


Analytic  Frequencies 

Galerkin  Approximates 

m  =  0 

m  =  0 

n  =  1 

266.05 

406.80 

603.83 

n  =  1 

266.06 

406.80 

603.83 

n  =  2 

531.53 

532.11 

924.29 

n  =  2 

531.53 

532.11 

924.29 

n  =  3 

541.03 

798.16 

1362.11 

n  =  Z 

541.03 

798.22 

1362.19 

n  =  A 

543.53 

1064.22 

1807.86 

II 

543.58 

1064.82 

1808.81 

n  =  b 

544.60 

1330.27 

2255.52 

n  =  5 

544.85 

1333.95 

2262.44 

m 

=  1 

m 

=  1 

II 

o 

338.75 

n  ==  0 

338.75 

n  —  \ 

147.09 

508.59 

873.61 

71=  1 

147.09 

508.59 

873.61 

n  =  2 

327.96 

714.55 

1115.55 

n  =  2 

327.96 

714.44 

1115.55 

n  —  Z 

434.81 

909.78 

1486.93 

n  =  3 

434.83 

909.85 

1487.00 

ra  =  4 

483.84 

1137.35 

1900.34 

II 

483.92 

1138.09 

1901.13 

n  =  5 

507.08 

1383.62 

2329.09 

n  =  5 

507.48 

1388.30 

2335.07 

m 

=  2 

m 

=  2 

O 

II 

677.50 

n  =  0 

677.50 

n  =  1 

64.56 

755.92 

1340.06 

n  =  1 

64.56 

755.92 

1340.06 

n  =  2 

187.34 

917.53 

1521.01 

n  =  2 

187.35 

917.53 

1521.01 

n  =  3 

296.26 

1100.21 

1804.82 

n  =  Z 

296.28 

1100.29 

1804.85 

n  =  4 

373.46 

1300.70 

2153.02 

n  =  4 

373.61 

1301.67 

2153.52 

n  =  5 

423.98 

1519.53 

2536.80 

n  =  5 

424.54 

1521.01 

2540.73 

m 

=  3 

m 

=  3 

n  =  0 

1016.25 

n  =  0 

1016.25 

n  =  1 

33.48 

1061.88 

1858.96 

n  =  1 

33.48 

1061.88 

1858.96 

ra  =  2 

111.09 

1177.72 

2001.04 

n  =  2 

111.09 

1177.73 

2001.04 

n  =  Z 

198.51 

1332.12 

2225.92 

n  —  Z 

198.54 

1332.21 

2225.94 

n  =  4 

275.56 

1509.99 

2514.74 

II 

275.76 

1511.08 

2515.01 

n  =  5 

336.44 

1706.31 

2848.26 

n  =  5 

337.29 

1714.42 

2850.19 

Table  5.  Frequencies  approximated  with  M  =  3  and  N  =  8  cubic  basis  functions. 
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m  =  0,n  =  5 


N  =  8  544.85  1333.95  2262.44 

iV  =  16  544.60  1330.30  2255.57 

Analytic  544.60  1330.27  2255.52 


m  =  l,n  =  5 


N  =  8  507.48  1388.30  2335.07 

N  =  16  507.08  1383.66  2329.14 

Analytic  507.08  1383.62  2329.09 


m  —  2,n  =  5 


II 

00 

424.54 

1521.01 

2540.73 

A  =  16 

423.99 

1519.56 

2536.84 

Analytic 

423.98 

1519.53 

2536.80 

m  =  3, 

71  =  5 

N  =  8 

337.27 

1714.42 

2850.19 

N  =  16 

336.45 

1706.35 

2848.28 

Analytic 

336.44 

1706.31 

2848.26 

Table  6.  Frequencies  approximated  with  M  =  Z  and  iV  =  8, 16  cubic  basis  functions 


U'^{x,Trj4) 


7r/4) 


4  e  B  10 


17^(a:,7r/4) 


V'^{x,Tr  /4) 


W^(x,Tr/4) 


Figure  4.  Approximate  873.61  Hz  modes  along  the  line  i  =  {(x,6)\6  =  •7r/4, 0  <  x  <  12} 


5.2  Modal  Example:  Fixed-Edge  Shell 

The  lack  of  analytic  expressions  for  the  frequencies  of  a  clamped-end  shell  necessitates  their 
numerical  approximation.  Like  the  simply  supported  case,  the  frequencies  are  approximated 
by  solving  the  matrix  eigenvalue  problem  (4.13)  with  basis  functions  constructed  so  as  to 
satisfy  the  fixed-edge  boundary  conditions  (2.6).  Axisymmetric  frequencies  obtained  using  a 
Fourier  limit  of  ilL  =  0  and  a:-axis  basis  number  of  =  16  are  summarized  in  Table  7.  A  com¬ 
parison  with  the  axisymmetric  frequencies  of  the  simply  supported  shell  (see  Table  3)  shows 
that  only  the  torsional  frequencies  match  (the  torsional  frequencies  are  boxed  to  facilitate 
comparison). 


Approximate  Frequencies 


266.05 

410.10 

532.11 

533.22 

541.76 

543.99 

544.96 

545.56 

546.07 

546.69 

547.63 

549.25 

552.18 

557.39 

565.57 

571.17 

575.20 

599.46 

798.16 

921.82 

1064.22 

1330.30 

1360.78 

1596.45 

1806.93 

1862.80 

2129.72 

2254.85 

2397.95 

2669.00 

2703.58 

2945.38 

3153.02 

3230.49 

3527.16 

3603.61 

3833.08 

4056.54 

4119.08 

4514.29 

4981.08 

5167.03 

5172.22 

5462.62 

5963.62 

6480.17 

8742.91 

6962.94 

8734.17 

Table  7.  Cubic  spline  approximates  of  axisymmetric  frequencies  obtained  with  M  =  0  and 
A  =  16  for  the  clamped-end  shell. 

While  the  torsional  modes  for  the  fixed-edge  shell  match  those  for  a  simply-supported 
shell  when  m  =  0  (see  [16]),  the  general  modes  for  the  two  shells  differ  due  to  boundary 
conditions.  To  illustrate,  fixed-edge  and  simply-supported  shell  modes  having  similar  fre¬ 
quencies  and  wavenumbers  are  compared  in  Figure  5.  Specifically,  the  545.56  Hz  mode  for 
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12  0  2  4  8 


10  12  0  2 


6  8  10  12 


545.66  Hz,  Simply  Supported 


Figure  5.  Axisymmetric  modes  U^{x,  0),  V^{x,  0)  and  W^{x,  0)  corresponding  to  the  fixed- 
edge  545.56  Hz  frequency  and  simply  supported  545.66  Hz  frequency. 


5.3  Forced  Shell  Examples 

A  second  means  of  testing  the  capabilities  of  the  approximation  method  is  through  the  consid¬ 
eration  of  steady  state  and  time-dependent  shell  responses  to  given  force  and  moment  inputs. 
In  the  first  examples  of  this  section,  known  true  solutions  are  used  to  calculate  corresponding 
force  input.  By  comparing  the  resulting  approximate  solutions  with  the  original  true  solu¬ 
tions,  the  expected  OiK^)  and  convergence  rates  for  the  linear  and  cubic  splines  can 

be  verified.  In  the  final  example,  a  triangular  input  to  the  patch  terms  (modeling  a  voltage  to 
spike  to  the  patches)  provides  a  broadband  model  response.  The  resulting  frequencies  are  then 
compared  with  those  obtained  by  solving  the  matrix  eigenvalue  problem  (4.13)  to  illustrate 
the  consistency  among  techniques. 
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5.3.1  Steady  State  Solution  —  Linear  and  Cubic  Spline  Approximations 

In  this  example,  steady  state,  ^-dependent  solutions  satisfying  the  clamped  edge  boundary 
condition  are  considered.  The  true  displacements  were  taken  to  be 

u{x,  6)  =  sin(27ra:/£)  sin(^) 

u(x,^)  =  sin(47rx/^)sin(20)  (5-1) 

w{x^0)  =  [cos(27ra:/.£)  —  1]  cos(30)  , 

and  hence  they  satisfy  the  fixed-edge  conditions  (2.6).  Plots  of  those  displacements  are  given  in 
Figure  6  to  illustrate  that  significant  bending  deformations  are  present.  The  forcing  functions 
were  obtained  by  substituting  the  true  solutions  into  the  strong  form  (2.1)  of  the  equations 
and  solving  for  the  functions  qe  and 

Since  the  forces  are  time-independent,  the  coefficients  =  [ui,  •  •  ■  ,ujv/'u5 
^1,  ■  •  ‘  for  the  approximate  displacements  u^(O,x),v^(0,x)  and  w^{9,x)  of  (4.1)  are 

obtained  by  solving  the  matrix  system 


=  F^  , 

where  and  are  defined  in  (4.9)  and  (4.10),  respectively.  The  specific  construction 
and  dimensions  of  and  F'^  depend  upon  the  bases  used  in  forming  and  (see 

Table  1). 

When  approximating  the  solutions  to  (5.1),  both  the  linear  and  cubic  axial  basis  functions 
summarized  in  Table  1  were  considered.  For  both  cases,  the  maximum  absolute  errors  over 
a  uniform  25  x  25  grid  in  (x,  0)  are  summarized  in  Table  8  while  ratios  of  the  errors  are 
given  in  Table  9  (the  absolute  errors  provide  the  same  information  as  relative  errors  since  the 
magnitudes  of  the  solutions  are  on  the  order  of  one).  For  both  sets  of  approximations,  equal 
axial  index  limits  were  assumed  and  are  denoted  by  —  N.  Due  to  the  nature 

of  the  forcing  function,  the  Fourier  limit  M  =  3  was  sufficient  for  resolving  the  circumferential 
behavior. 

The  error  results  in  Table  8  illustrate  the  rapid  convergence  of  the  method  with  the  fully 
cubic  spline  bcisis  while  the  ratios  in  Table  9  demonstrate  the  respective  O(h^)  and  0{h^) 
convergence  rates  of  the  cubic  and  linear  bases,  respectively  (the  multiplying  constant  causes 
the  ratios  in  the  linear  case  to  be  slightly  less  than  the  expected  value  of  four).  It  should  be 
noted  that  due  to  the  coupling  between  the  longitudinal,  circumferential  and  transverse  dis¬ 
placements  in  the  modeling  equations,  the  accuracy  of  the  transverse  approximate  is  reduced 
to  0{hl)  when  linear  splines  were  used  to  approximate  u  and  v,  even  though  cubic  splines 
were  used  to  approximate  w. 

The  approximate  longitudinal  displacements  obtained  with  N  —  4, 8, 16  and  32  linear 
axial  basis  functions  are  plotted  in  Figure  7.  Comparison  with  the  true  solution  in  Figure  6 
illustrates  that  with  N  =  16,  the  approximation  method  has  captured  the  qualitative  behavior 
of  the  solution.  Corresponding  plots  of  the  circumferential  and  transverse  displacements, 
plotted  in  Figures  14  and  15,  can  be  compared  with  the  true  circumferential  and  transverse 
solutions  in  Figure  6. 
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u 


V 


w 


Table  8.  Absolute  errors  in  the  longitudinal  u,  circumferential  u,  and  transverse  displace¬ 
ments  w  with  cubic  and  linear  spline  axial  basis  functions. 


Cubic  Basis 

Linear/ Cubic  Bcisis 

N=S 

N=4 

iV=8 

N=16 

N-8 

N=16 

N=S 

iV=16 

N=32 

u 

83.052 

50.888 

3.956 

3.863 

4.228 

V 

37.899 

24.069 

5.336 

3.563 

3.785 

w 

76.961 

73.395 

3.801 

3.966 

3.268 

Table  9.  Ratios  of  the  absolute  errors  with  linear  and  cubic  spline  axial  basis  functions. 


u{x,d)  v(x,9) 


Figure  6.  True  longitudinal,  circumferential  and  transverse  displacements. 


Approxlmat®  u,  N-4  Approximate  u,  N-8 


Figure  7.  Approximate  longitudinal  displacement  u  with  A'  =  4,  8, 16,32. 


Approximate  v,  N-4 


Apprmdmate  v,  N-S 


Approximate  v,  N-1 6 


Figure  8.  Approximate  circumferential  displacement  v  with  iV  =  4, 8, 16,32. 
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Approxlmai«  w,  N-4 


Approx)mal»w,N-8 


Figure  9.  Approximate  transverse  displacement  w  with  N  =  4,8, 16,32. 


5.3.2  Time-Dependent  Solution  —  Cubic  Spline  Approximation 

To  illustrate  the  approximation  of  the  dynamic  shell  model  (including  mass,  stiffness,  damping 
and  force  contributions),  time-dependent  forcing  functions  and  solutions  are  considered  in  the 
final  two  examples.  Cubic  spline  axial  basis  functions  are  used  for  discretizing  the  longitudinal, 
circumferential  and  tangential  displacements  in  both  cases. 

For  this  example,  the  true  solutions  were  taken  to  be 

u(x,  0)  =  sin(27ra:/^)  sin(^) 

v(x,  0)  =  sin(47rx/£)  sin(2^) 

w(x,0)  =  [cos(27rx/£)  —  1]  cos(30) 

(compare  with  (5.1)).  Solving  for  gx,ge  and  in  the  strong  form  of  the  equations  (2.1) 
yields  the  forcing  functions.  Since  the  forcing  terms  are  time-dependent,  discretization  of  the 
modeling  equations  yields  the  matrix  equation  (4.12)  in  R^,N  =  2A/’.  As  detailed  in  [16], 
the  resulting  ODE  system  is  moderately  stiff  and  a  variable-stepsize,  variable-order  stiff  ODE 
routine  was  used  to  numerically  integrate  the  system. 

The  absolute  errors  in  the  longitudinal,  circumferential  and  transverse  displacements  at 
time  T  =  1.5  sec  are  reported  in  Table  10  while  error  ratios  as  N  was  doubled  are  reported  in 
Table  11.  A  comparison  between  the  errors  in  Table  10  and  the  steady  state  errors  in  Table  8 
indicates  that  they  are  of  the  same  magnitude  while  comparison  between  Tables  11  and  9 
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indicates  nearly  identical  error  ratios.  Both  Table  10  and  11  comfirm  the  0{h‘^)  convergence 
rate  of  the  method  in  the  damped,  time-dependent  problem.  Hence  with  the  fully  cubic  spline 
axial  basis,  the  accuracy  contributes  to  the  efficiency  of  the  method  through  reduced  system 
sizes.  This  significantly  speeds  the  time  marching  which  is  advantageous  in  simulations, 
parameter  estimation  and  control  applications. 


Ar  =  4 

N  =  8 

N  =  ie 

u 

3.742  -  01 

4.506  -  03 

8.855  -  05 

V 

1.280  +  00 

3.376  -  02 

1.403  -03 

w 

7.056  -  01 

9.169  -  03 

1.249  -  04 

Table  10.  Absolute  errors  of  the  longitudinal  u,  circumferential  u,  and  transverse  displace¬ 
ments  w  at  time  T  =  1.5  seconds. 


oo 

II 

N  =  8 

TV  =  16 

u 

83.052 

50.888 

V 

37.899 

24.069 

w 

76.961 

73.400 

Table  11.  Ratios  of  the  absolute  errors  at  time  T  =  1.5  seconds. 


5.3.3  Time-Dependent  Solution  -  Patch  Input 

In  the  final  example,  the  external  inputs  to  the  patch  were  assumed  to  come  only  from  moments 
and  forces  generated  by  a  single  patch  pair  located  between  xi  =  4.5  in,  X2  =  7.5  in  and 
$1  =  0,  ^2  =  ’'■/b.  Material  and  patch  parameters  are  summarized  in  Table  12.  The  reported 
patch  parameters  are  consistent  with  those  of  patches  currently  being  used  in  the  Acoustics 
Division,  NASA  Langley  Research  Center. 

For  this  case,  the  force  vector  F^{t)  of  (4.10)  has  the  components 


r^2  r^2 

lF.it)h  =  l  I  R(N.Ut)^d.de 

= -  C  C 


dx  dO 


where  (iV^)pe,  (iVe)pe,  (Afa:)pe  and  (Me)pe  are  defined  in  (2.4)  with  i  =  1.  It  is  noted  that  due  to 
the  characteristic  functions,  the  integrals  reduce  to  only  that  region  covered  by  the  patches. 
This  is  in  accordance  with  the  physical  fact  that  moments  and  forces  are  generated  only  in 
patch  regions. 
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The  voltage  Vi(t)  to  the  outer  patch  (see  (2.5))  was  taken  to  have  a  triangular  shape 
with  a  magnitude  of  95  V  and  a  temporal  duration  of  1  x  lO""*  seconds.  The  voltage  V2{t)  to 
the  inner  patch  was  taken  to  be  V2{t)  =  This  input  models  an  out-of-phase  voltage 

spike  to  the  patches  and  serves  to  generate  a  broadband  shell  response.  The  input  of  such 
voltages  to  patches  is  a  current  technique  for  eliciting  shell  responses  when  performing  system 
identification  or  damage  detection. 

The  clamped-edge  boundary  condition  (2.6)  was  assumed  with  compatible  cubic  basis 
functions  in  the  numerical  approximation.  The  index  limits  N  =  8  and  M  =  3  were  sufficient 
for  resolving  the  shell  dynamics  in  the  previous  examples,  and  those  limits  were  also  used  to 
obtain  the  results  here.  The  frequencies  and  displacements  for  the  point  (x,  6)  =  (7.83, 0.89) 
are  plotted  in  Figures  10-12  (this  point  lies  just  outside  the  (x,0)  =  (7.5,7r/6)  corner  of  the 
patch).  Since  out-of-phase  simulated  voltages  were  input  to  the  inner  and  outer  patches,  the 
input  was  primarily  in  the  form  of  external  moments,  with  only  higher  order  terms  contributing 
to  the  in-plane  forces  (see  the  expressions  of  (2.5)).  The  small  magnitudes  of  the  displacements 
results  from  the  use  of  a  single  actuating  patch  pair  with  a  short  duration  voltage  spike.  Note 
that  all  computations  were  performed  with  double  precision  accuracy. 

In  Table  13,  the  natural  frequencies  obtained  through  this  dynamic  patch  excitation  are 
compared  with  those  obtained  for  the  undamped  shell  by  solving  the  matrix  eigenvalue  prob¬ 
lem  (4.13).  As  expected,  the  two  sets  of  frequencies  are  nearly  identical  since  the  same  dis¬ 
cretization  limits  were  used  in  both  cases.  The  slight  differences  are  due  to  the  Kelvin- Voigt 
damping  and  sample  rates  in  the  dynamic  simulation. 


h  =  0.01  in 
p  =  0.283  Ib/in^ 
z/  =  0.3 
i  =  12.0  in 

R  =  3.00  in 

F;  =  3x10^  Ib/in^ 

CD  =  15.099  lb  •  in  •  s 

X\  =  4.5  in 
=  0 

Eps  =  9.137  X 10®  Ib/in^ 
hpe  =  0.001  in 

X2  =  7.5  in 

02  =  7r/6 

Z/pe  =  0.31 

dsi  =  7.4x10-^  in/V 

Table  12.  Material  and  piezoceramic  patch  constants. 


t 
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Dynamically  Excited 
Frequencies  (Hz) 


Frequencies  from 
Eigenvalue  Problem  (Hz) 


60.1 

99.5 

60.2 

99.5 

176.1 

211.2 

176.0 

211.1 

301.8 

342.7 

301.8 

342.9 

343.0 

427.2 

342.9 

427.4 

435.8 

488.0 

435.8 

488.4 

508.4 

534.7 

508.7 

534.7 

545.0 

561.2 

544.9 

560.6 

758.1 

889.0 

758.2 

889.0 

919.5 

956.4 

919.3 

956.5 

119.4 

1119.2 

Table  13,  Natural  frequencies  obtained  through  dynamic  excitation  and  solution  of  the 
matrix  eigenvalue  problem  (4.13). 


Displacements  in  the  u  direction 


Excited  Frequencies  in  the  u  cfirection 
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Figure  10.  Longitudinal  displacements  and  frequencies  for  the  damped  shell. 


Displacements  in  the  v  direction 


Exdted  Frequencies  in  the  v  direction 


Figure  12.  Transverse  displacements  and  frequencies  for  the  damped  shell. 


6  Concluding  Remarks 

In  this  paper,  a  spHne-based  Galerkin  method  appropriate  for  discretizing  cylindrical  shell 
models  has  been  presented.  The  emphasis  throughout  the  discussion  was  on  the  development 
of  a  method  which  is  sufficiently  flexible  for  a  variety  of  models  including  those  which  arise 
in  smart  material  applications.  Galerkin  expansions,  in  terms  of  basis  functions  composed 
of  modified  splines  in  the  axial  direction  and  Fourier  components  circumferentially,  provided 
this  flexibility.  Discretization  in  this  manner  also  provided  an  accurate  and  efficient  means 
of  approximating  static  and  dynamic  shell  deformations  as  well  as  natural  frequencies  and 
modes  for  a  variety  of  boundary  conditions. 

As  demonstrated  through  numerical  examples,  the  method  yielded  0{h^)  convergence  for 
physically  realistic  values  of  the  shell  thickness  h  when  linear  splines  were  used  eis  axial  basis 
elements  for  u,v,  and  0{h^)  accuracy  when  all  three  components  u,v,w  were  approximated 
by  cubic  spline  expansions.  The  fourth  order  accuracy  attained  with  cubic  splines  contributes 
to  the  method’s  efficiency  by  reducing  the  number  of  basis  elements  and  the  hence  system 
sizes  necessary  for  resolving  the  shell  behavior  (in  many  of  the  examples,  N  =  S  cubic  spline 
basis  functions  were  sufficient  for  resolving  the  displacements  u,v,w).  Unlike  the  splines,  the 
Fourier  expansions  in  0  form  a  truncated  approximate  in  which  all  Fourier  components  with 
index  less  than  the  limit  index  M  are  resolved.  In  applications,  the  number  of  excited  Fourier 
modes  is  often  quite  small,  thus  permitting  small  limits  M. 

Another  advantage  of  the  cubic  spline  basis  is  the  fact  that  differentiability  and  smooth¬ 
ness  properties  are  constructed  in  the  basis  as  compared  with  many  finite  element  methods  in 
which  these  conditions  are  obtained  through  compatibility  conditions  at  the  nodes.  This  sim¬ 
plifies  construction  of  system  matrices  since  components  are  constructed  simply  through  the 
quadrature  of  basis  elements.  This  also  reduces  the  susceptibility  of  the  method  to  membrane 
and  shear  locking  since  there  are  no  elements  per  se  in  which  asymptotic  incompatibilities 
can  arise.  Alternatively,  the  spline-based  Galerkin  method  can  be  thought  of  as  similar  to  a 
spectral  method  but  with  basis  functions  having  finite  support. 
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Initial  tests  have  indicated  no  tendency  toward  membrane  or  shear  locking  when  using 
the  fully  cubic  spline  axial  bases  with  small  partition  size  h^.  This  is  consistent  with  finite 
element  results  reported  in  [25].  These  locking  considerations,  in  combination  with  the  0{h^) 
accuracy  of  the  cubic  splines,  indicate  that  the  fully  cubic  spline  axial  basis  will  be  preferable 
to  the  mixed  linear/cubic  basis  in  most  applications. 

The  method  was  also  demonstrated  to  be  flexible  with  regards  to  boundary  conditions, 
and  both  clamped-end  and  simply-supported  edge  conditions  were  attained  through  minor 
modifications  in  the  boundary  splines.  Free  end  conditions  are  natural  boundary  conditions 
in  the  weak  model  formulation  and  hence  no  basis  modifications  are  necessary  in  that  case. 
“Almost  fixed”  boundary  conditions  modeling  slight  edge  rotations  and  inplane  movement  are 
also  common  in  applications  due  to  the  difficulty  in  maintaining  perfect  boundary  clamps.  As 
discussed  in  [7],  only  minor  modifications  are  necessary  for  adopting  the  splines  for  this  case. 

Finally,  the  approximation  method  is  directly  applicable  to  shell  models  incorporating  pas¬ 
sive  (internal)  and  active  (external)  contributions  from  embedded  or  surface-mounted  sensors 
and/or  actuators.  When  derived  from  physical  principles,  the  passive  actuator /sensor  con¬ 
tributions  are  manifested  in  the  physical  model  parameters  (e.g.,  density,  stillness,  damping, 
Poisson  ratio)  and  are  incorporated  in  the  resulting  discrete  system  simply  through  piecewise 
integration  over  the  shell  domain.  The  external  moments  and  forces  are  isolated  to  regions 
covered  by  actuators  through  characteristic  functions  in  the  model.  In  the  weak  form,  this 
yields  integrals  over  the  actuator  regions.  While  partitions  should  be  aligned  with  actuator 
edges  to  obtain  optimal  accuracy,  neither  the  basis  nor  inner  products  must  be  modified  to  ac¬ 
count  for  sensors  and/or  actuators.  This  is  advantageous  when  approximating  shell  dynamics 
in  smart  material  applications. 
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